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1. Introduction 

Lattice QCD with staggered fermions in the strong coupling limit (SC-LQCD) is a useful effective 
model of QCD, as it shares important features of QCD such as confinement and spontaneous chiral 
symmetry breaking and its restoration at a transition temperature T c . One can study the nuclear potential 
as well as the phase diagram at non-zero baryon chemical potential Hb [1]. These topics can not be 
properly addressed with conventional, determinant-based lattice QCD using HMC algorithms, due to 
the notorious sign problem: all methods available today are limited to Hb/T < 1 [2]. In contrast, SC- 
LQCD can be reformulated as a monomer-dimer system [3]. There, the sign problem can be made mild 
due to aresummation of baryonic and mesonic degrees of freedom [4]. Due to algorithmic developments 
over the last decade, in particular due to the application of the Worm algorithm to the monomer-dimer 
partition function [5], SC-LQCD - which has been studied via mean field theory [6] and with Metropolis 
algorithms [3, 4] since the 1980s - has experienced a revival, as simulations at finite baryon density 
could be performed with modest computational demands. Moreover, the chiral limit can be studied very 
economically - simulations are faster than with a finite quark mass. However, limitations remain. In 
particular, only the 1 -flavor (4 tastes) theory has been considered so far in the dimer-formulation. The 
physically more interesting case of 2 flavors could not be addressed yet due to a severe sign problem 
in the mesonic sector [7]. Here, we propose a Hamiltonian formulation of strong coupling lattice QCD 
based on the Euclidean continuous time limit, where further simplifications occur. In particular, we show 
that this formulation is a generalization of Hamiltonians for spin systems. It can in principle be extended 
to arbitrary Nf. In this paper, we illustrate the formalism and give first Monte Carlo results obtained via 
Stochastic Series Expansion for Nf = 2 and U(2) gauge group. 



2. The continuous Euclidean time approach 



In SC-LQCD, the gauge coupling is sent to infinity and hence the coefficient j8 = 2N c /g 2 of the plaque- 
tte term representing the Yang Mills part F^ V F^ V of the action is zero. The lattice becomes maximally 
coarse, and no continuum limit can be considered. But the gauge fields in the covariant derivative can be 
integrated out analytically because the integration factorizes. After the Grassmann integration over the 
fermions, one obtains the SC-LQCD partition function [3] in the dimer representation, which is an exact 
rewriting of the 1 -flavor staggered fermion action on a d + 1 dimensional lattice N% x 7V T : 

S[U,x,X]=am q J £x(x)x{x) + l 1 £^) \x(x)e^U (x)x(x+6) -x(x+6) e - a ^tf (x)x(x) 



+ o EE \x(x)Ui(x)x(x+i)-x{x+t)uJ(x)x(x) 



(2.1) 



Z(m q ,n q) 



e n 



N c \k b \ 



{k,n,(} b=(x,H) 

w(£) = a(e)f^ eo{x)l exp(N c N t r e a T ix), (2.2) 

with m q the quark mass and jl = the quark chemical potential, a(£) = ±1 a geometry dependent 
sign and r{£) the winding number of baryon loop I. The sum Yl x% over admissible configurations, 
namely those which fulfill the Grassmann constraint 

N c 



n x + E 

fl=±d,...±d 



N c V* <E V. 



(2.3) 



Since color degrees of freedom have been integrated out, configurations are defined in terms of mesons 
- represented by the monomers n x G {0, . . .N c } and dimers (x) G {0, . . . N c } (non-oriented meson hop- 



2 



New algorithms and new results for SC-LQCD 



Wolfgang Unger 



pings) - and baryons - represented by self-avoiding closed loops constituted by £p(x) G {— 1,0, + 1}. 1 
Here, we consider the chiral limit, m q = where monomers are absent: n x = 0. 

In Eq. (2.1) we have introduced an anisotropy 7 in the Dirac couplings. This complication is necessary 
because the chiral restoration temperature is given by roughly aT ~ 1.5, and on an isotropic lattice 
with aT = 1 /N T we could not reach sufficiently high temperatures. Furthermore, varying 7 is the only 
way to vary the temperature continuously. The temperature is thus aT = f(y)/N r with f(y) = aja % . 
However, the functional dependence f(y) of the ratio of the spatial and temporal lattice spacings on 7 is 
not known. Naive inspection of the derivatives in Eq. (2.1) would indicate f(y) = 7, but this only holds 
at weak coupling. In contrast, the mean field approximation of SC-QCD based on a 1 /^-expansion [?] 
suggests that aT c = y 2 N t is the sensible, Af T -independent identification in leading order in l/d. We have 
emphasized elsewhere [8] by analytic arguments and numerical investigation that this identification is the 
only suitable one which renders observables like the chiral susceptibility and the specific heat finite in the 
limit N t , 7— > °°. However, this limit (keeping y 2 N t fixed) is approached with significant, sometimes non- 
monotonic 1 /N t corrections. To circumvent such extrapolation problems, we consider the continuous 
Euclidean time (CT) limit: N t — > °°, 7 — > °°, with y 2 /N x = aT fixed. Hence we are left with only one 
parameter /3 = N z /y 2 to set the thermal properties, and all discretization errors introduced by a finite Nx 
are removed. Moreover, in the baryonic sector the partition function simplifies greatly: baryons become 
static in the CT limit, hence the sign problem is completely absent. Additionally, multiple spatial dimers 
kj(x) > 1 become completely suppressed (see [8]) and one can derive the CT partition function: 

ZctGM)= K/3/2)* £ e^X L 4 T with k = £^ = ^4^, N L/T =Zn L/T (x) 

(2.4) 

where B is the baryon number, and is the set of equivalence classes W of graphs containing a total 
number k of spatial hoppings, equivalent up to time shifts of the vertices. The vertex weights vi = 1 
and vj = 2/\/3 label L- and T-types of vertices as illustrated Fig. 1 (left). An important property of the 
partition function Eq. (2.4) is that spatial dimers are distributed uniformly in time. The lengths A/3 of 
"dashed" or "solid" time intervals (see Fig. 1 left) are then, according to a Poisson process, exponentially 
distributed: P(Afi) oc exp(— AA/3), A/3 G [0,/3] with X the "decay constant" for spatial dimer emission 
X = (2d — Y,u ns(xdz A))/4. This is the basis for the continuous time Worm algorithm presented in [8]. 

3. The 1-flavor Hamiltonian 

The Hamiltonian formulation can be obtained from Eq. (2.4) by realizing that the degrees of freedom 
can be mapped on a spin system. We can restrict the discussion to the mesonic sector U(N C ), since 
baryons are static for N c > 2. Notice that, except for spatial hoppings, meson lines are time-like and 
form dimer chains alternating between k^ en (x) G {0, . . .N c } and k^ d (x) with k% dd (x) = (N c - ^ ven (x)) 
on even and odd time-slices. One can then introduce the observable 

S z {x,t) = y - '— (2ko(x,t)-N c ) G {-N c /2,...N c /2} (3.1) 

which is constant on static lines. The "spin" S z simply counts the number of time-like meson hoppings, 
and is in no way related to the spin of the quarks. Spatial dimers can then be oriented consistently, 
such that for each spatial dimer between a pair of neighboring sites (x,y), one unit of spin AS^ = ±1 is 
transferred from site x to site y. Hence the total spin S z = S z (t) = Y,xeV S z (x,t) is globally conserved. 

'Note that U(3) describes a purely mesonic system, while SU(3) contains baryons. 
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Figure 1 : Left: Typical 2-dimensional configurations in discrete (top) 
and continuous time (bottom) at the same temperature, where multiple 
dimers are absent and baryons (red) become static. L- and T-vertices 
have different weights. Meson lines can be oriented in a consistent 
way. Right: Mapping of static mesons (time-like chains) to a "spin". 



In the case of gauge group U(l) (lattice QED in the strong coupling limit), the Hamiltonian is identical 
to that of the XY Model in zero field: H = L{x.v) G + G > with = G\ ± ia% the spin raising/lowering 
operators constructed from the Pauli matrices. The generalization to U (N c ) reads 



Z(j3)=Tr 



H 



1 



Jp^+h.c. 



vi 



V 



Vk 



lk(l+N c -k) 



vn c 0, 



(3.2) 

with J = (J + ) T which are spin lowering/raising operators. The off-diagonal matrix elements Vk are 
generalized vertex weights. For Af c = 3, we can identify vz, = vi = V3 = 1, Vj = V2 = 2/y/3. Note that 
the operators reflect the existence of a lowest and highest weight, J \ —N c /2) = 0, J + \N c /2) = 0, and 
fulfill as well the commutation relation ^[/ + ,/~] = J z = diag (-N c /2, . . .N c /2) with J Z \S Z ) =S Z \S Z ). 
This justifies the characterization of static lines in terms of a "spin" quantum number. 



4. Stochastic Series Expansion 

The partition function Eq. (2.4), which is a sum of weighted diagrams in a perturbative series in j8, 
can be sampled via diagrammatic Monte Carlo techniques such as the continuous time worm algorithm 
[9] or the Stochastic Series Expansion (SSE) [10]. Here, we restrict to SSE: as we will see, it can be 
easily generalized to Nf > 1 once we have constructed the corresponding Hamiltonian. SSE is based on 
a rewriting of the partition function by inserting identity and diagonal matrix elements: 



Z(j3)=Tr 



,-PH 



II 

X {S L } 



P k (L-k)\ 
L! 



X 



1=1 



X 



Hih — ct, H 2 .b 



2 ^' 



(4.1) 

where % is a state vector and St is a time-ordered sequence of indices: Sl = {[ai,b\], [02,^2]) ■ • • [ol^l]} 
characterizing - together with an initial state % - a graph in Z(j8). L is the number of operators and 
K < L the order in the expansion in j8. The indices a, = correspond to the identity, a, = 1 to the 
diagonal matrix element cl, where c can be adjusted in order to simplify the algorithm, and a, = 2 to 
non-diagonal matrix elements. The index bj = (x,y) G Vd denotes a bond on the lattice, and for a, = 0, 
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bi = denotes a dummy bond. For any finite volume and given temperature, only a finite number of 
orders in j3 contribute. L can be set to be larger than this number, making SSE approximation-free. The 
algorithm consists of two kinds of updates: (1) a Metropolis update changing the order in j8: 

Nldp(x\ft hb \x) _ L-K+l 



P([0,0]-> [l,b]) 



P([l,b] ^ [0,0]) 



(4.2) 



and (2) the operator loop update, visiting a set of bonds in succession, starting from an input leg and de- 
termining the output leg with heatbath probability oc ^x{x)x{y)\H aj b j \x' (x)x' (y)) ■ In Fig. 2 we compare 
different observables obtained from SSE and/or CT-Worm. We want to stress that a new observable, the 
spin susceptibility Xs= y ((Li ^) 2 ), is also sensitive to the chiral transition. It measures the fluctuations 
in the number of time-like mesons, and is thus analogous (but not equal) to the specific heat. 
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Figure 2: Observables in 1-flavor U(3): (top row) energy density and susceptibility measured with CT- 
worm and SSE agree; (bottom left) spin susceptibility measured with SSE; (bottom right) chiral suscep- 
tibility measured with CT-Worm. 



5. The 2-flavor formulation and first results for U(2) gauge group 

In the strong coupling limit, only one of the four tastes of staggered fermions remains light, since taste 
splitting is maximal. It was argued in [1] that the nuclear interaction in 1-flavor SC-LQCD is due to 
entropic forces: the presence of static baryons modifies the pion bath, with the modulation proportional 
to the p -propagator. However, with 1 flavor pion exchange cannot not occur, because mesons do not 
couple to baryons. Pion exchange can only be studied by going to Nf > 1. But no 2-flavor formulation 
of staggered SC-QCD suitable for finite-density Monte Carlo exists, as already the mesonic sector has 
a severe sign problem [7]. After reviewing the problems with the conventional dimer representation of 
SC-QCD, we derive a sign-problem-free Hamiltonian formulation in continuous Euclidean time. 

The 1-link integrals which appear in the strong coupling limit can be expressed in terms of the gauge- 
invariant terms [11]. For SU(2) and U(2), the 1-link integals (valid for any Nf) are 

f (X + A)» f X'+ZJ 

= I^TI)!' M 2 )(^y)=£ o( . +2 . +1)! , ( . !)2 (3-D 
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Figure 3: Left: definition of unflavored, flavored, and mixed dimers. Right: 2x2 sample configurations 
with positive and negative weight. 



with X = tr(mrrfi) = U +D+7i + + 7T~ the sum of Nf = 2 unflavored (\J=uu andD=Jfi?) andNf(Nf— 1) =2 
flavored mesons, A = det(m) + det(m t ) the diquark+anti-diquark term, and Z = det(mra^) = X 2 — UD — 
7t + 7t~ + a + j6 a mesonic term. Here, a = ud x du x iiu y dd y and j8 = uu x dd x ud y du y are the potentially 
problematic contributions, which correspond to the mixing of UD and 7T + 7r~ dimer pairs (obtained via 
non-trivial Wick-contractions), as illustrated Fig. 3. In particular, if a configuration contains an odd 
number of a or j3 links, according to the Grassmann constraint (see [7]) the configuration has a negative 
sign. The essential feature of the continuous Euclidean time formulation is that multiple spatial dimers, 
and hence also a and j3 spatial dimers, are suppressed. They can only enter in static lines, where they 
can be resummed with other static lines so that the sign problem disappears completely. Combining 
time-like dimers of alternating orders is analogous to the procedure discussed in the Nf = 1 case [8]: 
Chains of alternating orders &(y lk ) x 0[y w ^ N i^ lk ^ are resummed in a way consistent with the constraint 
Eq. (2.3). This gives rise to new conserved quantum numbers: the spin is now composed of iVf separate 
spins S z = Z%i Sf e {-AWc/2, . . .NfN c /2}, with 5^ G {-N c /2, . . .N c /2} each being conserved, and 
we also get Nf(Nf — 1)/2 charges <2,- G {—N c ,. . . , +-/V C }- In the specific case of two flavors (/ = U,D and 
<2i = <2rc)> tne spins Su,Sd can be replaced by S z and Qi = Sfj — Sj) G {— N c , +N C } which may be 
viewed as isospin. Here, we find in total 19 types of states, as illustrated in Fig. 4. The state vector and 
the transition rules at spatial dimers are given by 



X 



{5^(x)}/=i,...iV f ,{!2i(^)}/=l,...iV f (JVf-l)/2) , 



\AS Z 



1 and £|Afi(|€{0,l}. (5.2) 



The Hamiltonian for Nf = 2 is now a sum of four contributions, implementing these transition rules: 
1 



H - ? XI \ J U(x) J U(v) +J D(x) J D(y) +J n 



+ {x) J 7t+(y) +J M-(x) J M-(y) +hX - 



(5.3) 
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Figure 4: Classification of all 
static lines for U(2) based on 
"spin" quantum number and 
particle content. Each static 
line consists of two arrows, 
dashed arrows denoting spin, 
curly arrows denoting flavor- 
neutral mesons content, and 
solid arrows denoting charged 
mesons. 



2 For SU(2) we also have 3 kinds of diquarks, uu, dd and ud, which are not suppressed in the CT-limit. To avoid this 
complication, for N c = 2 we restrict to U(2). For N c = 3 this restriction is not necessary. 
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The absorption /+, and emission 7^ operators can be represented as 19 x 19 lower left/upper right 
triangular matrices, where the entries are again given by vertex weights: = ^= if states with S = 
and \Q\ = 1 are connected to the vertex, v Kj = 1 otherwise. This Hamiltonian can be used in the SSE 
algorithm to obtain results on the chiral phase transition. Our preliminary Nf = 2 results are compared 
to mean field predictions in Tab. 1 . 



N c Nf = 1 Nf = 2 Table 1 : Comparison of the critical temperature aT c between 

1 3/2 [1.102(1)] 5/5 [0.77(1)] mean field results and Monte Carlo results [in brackets] for 

2 4/2 [1.467(1)] 6/5 [1.04(1)] U(7Y C ) gauge groups. The new results are in column Nf = 2. 

3 5/2 [1.884(1)] 7/5 The MC value for Nf = 2, U(3) has not been measured yet. 



6. Conclusion 

We have given a new, Hamiltonian formulation of strong coupling lattice QCD with staggered fermions 
in the chiral limit. It is based on the insight that strong coupling lattice QCD in the continuous time limit 
is analogous to a spin system. A new observable, the spin susceptibility, turns out to be sensitive to the 
chiral transition. Also, the Hamiltonian description allows to apply quantum Monte Carlo methods. In 
[8] we have studied 1-flavor thermodynamics via the continuous time Worm algorithm. Here, we make 
the first step towards 2-flavor simulations, by making use of Stochastic Series Expansion, a diagrammatic 
Monte Carlo technique which we generalize to the Hamiltonian in question. SSE has the advantage that 
more complicated Hamiltonians can be studied with ease. The drawback of SSE, in contrast to the 
continuous time Worm, is that we do not know of a way to obtain the two-point function for free. For the 
computation of the specific heat, the performance of both algorithms is quite similar. We have provided 
first results on the U(2) transition with two flavors. The extension to SU(3) with finite baryon chemical 
potential is more involved: the number of static lines increases to 44 in the mesonic sector. 
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